TP 1 : Estimation
🍰 Estimation du nombre de boulangeries dans le Nord.
Total
Un statisticien souhaite déterminer le nombre de boulangeries dans le département du Nord. Pour faire cela, il souhaite effectuer un sondage dans lequel il tirera des communes parmi les 638 appartenant au département et demandera le nombre total de boulangeries dans la commune.
Dans un premier temps, le statisticien décide d’utiliser un plan aléatoire simple sans remise de taille 60.
Les résultats sont disponibles dans le fichier suivant : ech_srs.csv.
- Décrivez la population, la variable d’intérêt et la fonction d’intérêt. Dans la suite, nous noterons \(\mathcal{U}\), la population et \(y_k\) la valeur de la variable d’intérêt pour l’individu \(k\). > La population est l’ensemble des communes du Nord, la variable d’intérêt (pour une commune) est le nombre de boulangeries au sein de la commune et la fonction d’intérêt est le nombre total de boulangeries dans la commune.
- Proposez un estimateur sans biais du total du nombre de boulangeries dans le département du Nord. > L’échantillon est tiré selon un sondage aléatoire simple sans remise de taille 60 : les probabilités d’inclusion d’ordre 1 sont non nulles pour tout individu de la population. Il en vient que l’estimateur du total d’Horvitz-Thompson est un estimateur sans biais du total.
- Proposez une estimation associée à cet estimateur (on pourra utiliser la fonction
HTestimation).
library("sampling")
ech <- read.csv("ech_srs.csv",
sep = ",",
colClasses = c("COG" = "character"))
#Vérification de la table
str(ech)
n <- 60
N <- 638
#Calcul des probabilités d'inclusion d'ordre 1
pik_srs <- rep(n/N, n)
tot_boul_srs <- HTestimator(y = ech$BOULANGERIES, pik = pik_srs)
cat("L'estimation du total basée sur l'estimateur d'Horvitz-Thompson avec un échantillon tiré selon un SRS donne : ", tot_boul_srs, " boulangeries.")- Proposez un estimateur de la variance de l’estimateur proposé à la question 2 (on pourra utiliser la fonction
varHT).
#Calcul des probabilités d'inclusion d'ordre 2
pi_kl <- matrix(n*(n-1)/(N*(N-1)), ncol = n, nrow = n)
# Par contre il faut remplacer les termes diagonaux car pi_kk = pi_k
diag(pikl) <- n/N
#Methode 1 : estimateur de Horvitz-Thompson de la variance de l'estimateur du total de HT
estim_var_tot_boul_srs <- varHT(y = ech$BOULANGERIES, pi_kl, method = "1")
#Methode 2 : estimateur de Sen-Yates-Grundy de la variance de l'estimateur du total de HT
varHT(y = ech$BOULANGERIES, pi_kl, method = "2")
#Même résultat car SRS.- Donnez un intervalle de confiance asymptotique au niveau 0.90 du total. Calculez une réalisation de cet intervalle. Commentez.
#I
niveau <- 0.90
alpha <- 1 - niveau
quantile_norm <- qnorm(p = 1- (alpha/2))
#Borne inférieure :
tot_boul_srs - quantile_norm * sqrt(estim_var_tot_boul_srs)
#Borne supérieure :
tot_boul_srs + quantile_norm * sqrt(estim_var_tot_boul_srs)Le statisticien se rend compte tardivement qu’il peut se baser sur des informations complémentaires pour choisir les probabilités d’inclusion d’ordre 1 de son plan de sondage. En se rendant sur le (merveilleux) site de l’Insee, il trouve le nombre d’habitants pour chaque commune du Nord.
Il décide cette fois d’opter pour un plan poissonien en attribuant une probabilité d’inclusion d’ordre 1 plus grande pour les communes les plus peuplées. Les probabilités d’inclusion d’ordre un sur l’ensemble de la population sont comprises entre 0.00105398 et 1.
Les résultats sont disponibles dans le fichier suivant : ech_poisson.csv
- Quelle est la taille de l’échantillon ?
ech_poisson <- read.csv("ech_poisson.csv", sep = ",", colClasses = c("COG" = "character"))
nrow(ech_poisson)- Proposez un estimateur sans biais du total du nombre de boulangeries dans le département du Nord.
library(dplyr)
ech_poisson |>
filter(pi_k == 0) |>
count()- Proposez une estimation associée à cet estimateur.
tot_boul_poisson <- HTestimator(y = ech_poisson$BOULANGERIES, pik = ech_poisson$pi_k)- Proposez un estimateur de la variance de l’estimateur proposé à la question 2.
pi_kl_poisson <- t(t(ech_poisson$pi_k)) %*% ech_poisson$pi_k
diag(pi_kl_poisson) <- ech_poisson$pi_k
estim_var_tot_boul_poisson <- varHT(ech_poisson$BOULANGERIES, pikl = pi_kl_poisson, method = "1")- Donnez un intervalle de confiance asymptotique au niveau 0.90 du total. Calculez une réalisation de cet intervalle. Commentez.
tot_boul_poisson - quantile_norm * sqrt(estim_var_tot_boul_poisson)
tot_boul_poisson + quantile_norm * sqrt(estim_var_tot_boul_poisson)- Pourquoi de telles différences avec le plan de sondage de la première partie de l’exercice ? On pourra utiliser l’intuition fournie par la formule de Sen-Yates-Grundy (même s’il ne s’agit pas d’un plan de taille fixe dans cette partie).
Comme indiqué, l’information auxiliaire dont nous disposons pour chaque commune est le nombre d’habitants.
Cette variable étant positive, il est possible de calculer des probabilités proportionnelles au nombre d’habitants.
Une contrainte imposée ici est que \(\sum_{k \in \mathcal{U}} \pi_k = 60\) : autrement dit, l’espérance de la taille de l’échantillon sous le plan poissonien est égale à celle de l’échantillon tiré selon le SRS de la première partie de l’exercice.
Intuitivement, si \(\text{pop}_k\) désigne la population de la commune \(k\), on peut souhaiter que \(\pi_k = n \frac{\text{pop}_k}{\sum_{j \in s} \text{pop}_j}\). Néanmoins, certaines unités peuvent avoir une probabilité d’inclusion d’ordre 1 supérieure à 1. Dans ce cas, on applique l’algorithme suivant :
- \(A = \emptyset\)
- Calculez pour tout \(k \in \mathcal{U}/A\), \(\pi_k = (n - (Card(A))) \frac{\text{pop}_k}{\sum_{j \in s} \text{pop}_j}\).
- Soit \(b = \{k \in \mathcal{U}/A | \pi_k > 1 \}\) \(A \gets A \cup b\) Revenir à l’étape 2 tant que \(Card(b) \neq 0\).
Nombre moyen de boulangeries.
Le statiscien veut également fournir une information sur le nombre de boulangeries moyen \(\mu_y := \frac{1}{N} t_y\) dans le Nord en utilisant les données déjà récoltées.
Il opte pour deux stratégies :
- utilisation simple principe de subtitution : un estimateur \(\hat{\mu}_y^{(1)}\) est fourni en remplaçant \(t_y\) par l’estimateur d’Horvitz-Thompson du total \(\hat{t}_y\).
- utilisation double du principe de subtitution : un estimateur \(\hat{\mu}_y^{(2)}\) est fourni en remplaçant \(t_y\) et \(N\) (même s’il est connu) par des estimateurs sans biais.
Nous considérons ici uniquement les données issues du plan poissonien.
- Montrez que \(\hat{\mu}_y^{(1)}\) est un estimateur sans biais pour \(\hat{\mu}_y\) et proposez un estimateur de la variance de \(\hat{\mu}_y^{(1)}\)
- En notant que \(\displaystyle N = \sum_{k \in \mathcal{U}} z_k\) où pour tout \(k \in \mathcal{U} ~~ z_k = 1\), proposez un estimateur \(\hat N\) sans biais de N.
- Pourquoi il n’est pas possible de calculer explicitement un estimateur de la variance de \(\hat{\mu}_y^{(2)}\) ?